surge analysis on IOC tide gauge¶
case study: cyclone FUNG WONG
Install¶
Libraries needed for this notebook:
pip install searvey hvplot utide ipykernel geoviews pyarrow
In [1]:
import searvey
import utide
import pandas as pd
import hvplot.pandas
from shapely.geometry import box
/home/tomsail/work/gist/ioc_surge_analysis/.venv/lib/python3.12/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html from .autonotebook import tqdm as notebook_tqdm
In [2]:
ioc_stations = searvey.get_ioc_stations()
In [3]:
bbox = box(100.0, -5.0, 130.0, 25.0) # lon_min, lat_min, lon_max, lat_max
In [4]:
selected_stations = ioc_stations[ioc_stations.within(bbox)]
selected_stations
Out[4]:
| ioc_code | gloss_id | country | location | connection | added_to_system | observations_arrived_per_week | observations_expected_per_week | observations_ratio_per_week | observations_arrived_per_month | ... | transmit_interval | lat | lon | contacts | dcp_id | last_observation_level | last_observation_time | delay | interval | geometry | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 45 | ambon | 68 | Indonesia | Ambon | SZXX01 | 2008-11-20 12:43:00 | 9180 | 10080.0 | 91 | 36510 | ... | 5' | -3.683 | 128.183 | Geospacial Agency of Indonesia ( Indonesia ) +... | 301434060409970 | 5.03 | 13:00 | 9' | 5' | POINT (128.183 -3.683) |
| 47 | AMTSI | <NA> | Indonesia | Port of Ampana, Central Sulawesi | web | 2025-08-15 14:11:20 | 9690 | 10080.0 | 96 | 37857 | ... | 5' | -0.847 | 121.601 | Kepala Badan Meteorologi, Klimatologi, dan Geo... | NaN | 0.01 | 13:00 | 9' | 5' | POINT (121.601 -0.847) |
| 132 | BETSI | <NA> | Indonesia | Port of Beo, North Sulawesi | web | 2025-08-15 14:26:13 | 9835 | 10080.0 | 98 | 43470 | ... | 5' | 4.230 | 126.788 | Kepala Badan Meteorologi, Klimatologi, dan Geo... | NaN | 0.01 | 13:05 | 4' | 5' | POINT (126.788 4.23) |
| 136 | bitu | 69 | Indonesia | Bitung | SZXX01 | 2008-11-20 13:48:00 | 9125 | 10080.0 | 91 | 36675 | ... | 5' | 1.439 | 125.190 | Geospacial Agency of Indonesia ( Indonesia ) +... | 300434067058990 | 8.48 | 13:00 | 9' | 5' | POINT (125.19 1.439) |
| 183 | bupo | <NA> | Indonesia | Bungus Port | bgan | 2020-07-07 09:10:02 | -down- | 10080.0 | 0 | -down- | ... | 1' | -1.030 | 100.396 | Joint Research Centre ( Europe ) +Ministry of ... | BUPO | 1.45 | -down- | 1177d | 1' | POINT (100.396 -1.03) |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 1621 | txil | <NA> | Taiwan | Xiao Liuqiu | web | 2024-01-02 15:35:20 | 168 | 168.0 | 100 | 743 | ... | 1h | 22.353 | 120.383 | Central Weather Administration ( Taiwan ) | NaN | 0.41 | 12:00 | 1h | 1h | POINT (120.383 22.353) |
| 1622 | tyon | <NA> | Taiwan | Yongan | web | 2024-01-02 15:35:20 | 168 | 168.0 | 100 | 743 | ... | 1h | 22.819 | 120.198 | Central Weather Administration ( Taiwan ) | NaN | -999.00 | 12:00 | 1h | 1h | POINT (120.198 22.819) |
| 1672 | vung | <NA> | Viet Nam | Vung Tau | SZXX01 | 2007-10-27 11:39:00 | 8895 | 10080.0 | 88 | 35835 | ... | 6' | 10.340 | 107.071 | Hydro-meteorological and Environmental Station... | 300434067059990 | 3.81 | 13:00 | 9' | 6' | POINT (107.071 10.34) |
| 1706 | xish | <NA> | China | Xisha | SZCI01 | 2014-10-02 13:42:13 | -down- | 10080.0 | 0 | -down- | ... | 5' | 16.830 | 112.330 | China Meteorological Administration ( China ) | 11745 | 1.30 | -down- | 3134d | 5' | POINT (112.33 16.83) |
| 1717 | zhap | 78 | China | Zhapo | SZCI01 | 2014-10-02 14:12:15 | 4500 | 10080.0 | 45 | 18491 | ... | 1' | 21.580 | 111.820 | China Meteorological Administration ( China ) | 09731 | 2.59 | 12:59 | 10' | 1' | POINT (111.82 21.58) |
70 rows × 25 columns
In [5]:
selected_stations.hvplot.points(
geo=True,
tiles=True,
hover_cols=["ioc_code"],
title = "stations in the selected region"
)
Out[5]:
detide the stations¶
In [6]:
def surge(ts: pd.Series, lat: float, rsmp: int = None):
ts0 = ts.copy()
OPTS = {
"constit": "auto",
"method": "ols", # ols is faster and good for missing data (Ponchaut et al., 2001)
"order_constit": "frequency",
"Rayleigh_min": 0.97,
"lat": lat,
"verbose": True,
}
if rsmp is not None:
ts = ts.resample(f"{rsmp}min").mean()
ts = ts.shift(freq=f"{rsmp / 2}min")
coef = utide.solve(ts.index, ts, **OPTS)
tidal = utide.reconstruct(ts0.index, coef, verbose = OPTS['verbose'])
return pd.Series(data=ts0.values - tidal.h, index=ts0.index)
download the data¶
create a data folder with raw and surge subfolders
In [7]:
! mkdir -p data/{raw,surge}
let's first download all the station in the bounding box
In [13]:
drop_columns = ["bat"]
def serialize(d): # to export metadata
out = {}
for k, v in d.items():
if v is not pd.NA and not pd.isna(v):
out[k] = str(v)
return out
In [19]:
for irow, row in selected_stations.iterrows():
lat = row.lat
ioc_code = row.ioc_code
df_raw = searvey.fetch_ioc_station(
ioc_code,
pd.Timestamp.now()-pd.Timedelta(days=90), # we need at least 90 days (in theory we'd need more..) to remove properly tidal constituents
pd.Timestamp.now()
)
df_raw.attrs = {**serialize(dict(row)), **{"signal_type": "raw"}}
df_raw.to_parquet(f"data/raw/{ioc_code}.parquet")
IOC-bupo: No data. Creating a dummy dataframe IOC-ms002: No data. Creating a dummy dataframe IOC-ms004: No data. Creating a dummy dataframe IOC-ms005: No data. Creating a dummy dataframe IOC-ms006: No data. Creating a dummy dataframe IOC-nans: No data. Creating a dummy dataframe IOC-padn: No data. Creating a dummy dataframe IOC-sebla: No data. Creating a dummy dataframe IOC-tanjo: No data. Creating a dummy dataframe IOC-tdaw: No data. Creating a dummy dataframe IOC-tfan: No data. Creating a dummy dataframe IOC-thep: No data. Creating a dummy dataframe IOC-tluk: No data. Creating a dummy dataframe IOC-tnan: No data. Creating a dummy dataframe IOC-xish: No data. Creating a dummy dataframe
In [20]:
for irow, row in selected_stations.iterrows():
df_raw = pd.read_parquet(f"data/raw/{row.ioc_code}.parquet")
df_raw.loc["2025-11-01":].dropna().hvplot(
title = f"detided signal for {row.ioc_code}, {row.location}, lat:{row.lat}, lon:{row.lon}"
)
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
Out[20]:
bokeh backend could not plot any Elements in the Overlay.
Out[20]:
:NdOverlay [Variable]
Out[20]:
detide the station¶
Let's look at the stations of interest
In [21]:
detide_selection = {
"lega" :"rad",
"luba" :"prs",
"quin" :"ras",
"qing" :"rad",
"quar" :"flt",
"mani" :"prs",
"subi" :"rad",
"thsi" :"rad",
"txil" :"rad",}
In [22]:
for irow, row in selected_stations.iterrows():
if row.ioc_code in list(detide_selection.keys()):
df_raw = pd.read_parquet(f"data/raw/{row.ioc_code}.parquet")
df_surge = surge(df_raw[detide_selection[row.ioc_code]].dropna(), lat=row.lat, rsmp=2)
df_surge.loc["2025-11-01":].hvplot(
title = f"detided signal for {row.ioc_code}, {row.location}, lat:{row.lat}, lon:{row.lon}"
)
df_surge.attrs = dict(row)
df_surge.attrs = {**serialize(dict(row)), **{"signal_type": "detide"}}
df_surge.to_frame().to_parquet(f"data/raw/{row.ioc_code}.parquet")
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
solve: matrix prep ... solution ... done. prep/calcs ... done.
Out[22]:
In [23]:
ioc_stations[ioc_stations.ioc_code.isin(list(detide_selection.keys()))].hvplot.points(
geo=True,
tiles=True,
c = 'r',
s=100,
hover_cols=["ioc_code"],
title = "stations that recorded a surge"
)
Out[23]:
look at a specific station¶
In [24]:
station = "lega"
for irow, row in selected_stations.iterrows():
if row.ioc_code == station:
df_raw = pd.read_parquet(f"data/raw/{row.ioc_code}.parquet")
df_raw
station
df_raw.loc["2025--01":].dropna().hvplot(
title = f"signal for {row.ioc_code}, {row.location}, lat:{row.lat}, lon:{row.lon}"
)
Out[24]:
| 0 | |
|---|---|
| time | |
| 2025-08-13 14:16:00 | 0.172636 |
| 2025-08-13 14:17:00 | 0.172863 |
| 2025-08-13 14:18:00 | 0.174116 |
| 2025-08-13 14:19:00 | 0.173393 |
| 2025-08-13 14:20:00 | 0.174694 |
| ... | ... |
| 2025-11-11 13:01:00 | -0.149741 |
| 2025-11-11 13:02:00 | -0.149008 |
| 2025-11-11 13:03:00 | -0.148250 |
| 2025-11-11 13:04:00 | -0.148467 |
| 2025-11-11 13:05:00 | -0.147660 |
106149 rows × 1 columns
Out[24]:
'lega'
Out[24]: